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INTRODUCTION 

The  promising  potential  of  MRI  in  diagnosis  of  breast  cancer,  as  a  complementary 
modality  to  X-ray  mammography,  has  been  well  recognized  [1][2][3].  Contrast- 
enhanced  (CE)  MRI  is  regarded  as  the  most  sensitive  technique  for  breast  cancer 
detection;  in  addition,  it  is  extraordinarily  useful  in  assessing  the  tumor  extent  and 
assessing  the  therapy  response.  Despite  its  well-recognized  utilities,  however,  the 
technique  has  not  been  introduced  to  routine  clinical  breast  imaging.  One  of  the  most 
important  obstacles  is  the  lack  of  standardization  in  terms  of  interpretation  guidelines 
[2] [4],  The  reproducibility,  effectiveness  and  relative  significance  of  interpretation 
criteria  in  the  literature  are  far  from  being  well  evaluated.  The  purpose  of  the  proposed 
research  is  to  develop  computerized  methods  to  take  full  advantage  of  the  wealth 
information  that  dynamic  MRI  offers  to  improve  methods  for  the  diagnosis  and  prognosis 
of  breast  cancer.  The  research  involves  investigation  of  automatic  methods  for  image 
artifacts  correction,  tumor  segmentation,  and  extraction  of  computerized  features  that 
help  distinguish  between  benign  and  malignant  lesions.  Our  hypothesis  is  that  advanced 
methods  for  the  computerized  interpretation  of  state-of-the-art  breast  MRI  will  grant  MRI 
a  more  leading  role  in  diagnosis  and  monitoring  of  treatment  for  breast  cancer,  and 
greatly  benefit  breast  cancer  patients. 
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BODY 

Training  Accomplishments 

At  the  time  of  this  report,  the  recipient,  Weijie  Chen,  of  the  Predoctoral  Traineeship 
Award  has  taken  all  the  22  required  courses  towards  the  Ph.D.  degree  in  medical  physics. 
The  courses  include  physics  of  medical  imaging,  physics  of  radiation  therapy, 
mathematics  for  medical  physicists,  image  processing,  statistics,  machine  learning, 
numerical  computation,  topics  in  computer  vision,  anatomy  of  the  body,  radiation  biology, 
beam  physics,  clinical  physics  in  PET,  and  teaching  assistant  training. 

The  recipient  has  also  been  working  extensively  on  his  Ph.D.  dissertation  proposal, 
which  is  a  critical  step  towards  the  Ph.D.  degree  as  required  by  the  Committee  on 
Medical  Physics,  the  University  of  Chicago. 

Research  Accomplishments 

1.  Database  collection 

The  first  part  of  our  work  has  been  collecting  breast  MR  data.  We  have  previously 
collected  a  database  (referred  to  as  Database-I  in  the  following)  containing  images  from 
121  patients.  Images  were  obtained  using  a  T1 -weighted  3D  spoiled  gradient  echo 
sequence.  After  the  acquisition  of  the  precontrast  series,  Gd-DTPA  contrast  agent  was 
injected  intravenously  by  power  injection  with  a  dose  of  0.2mmol/kg.  Five  postcontrast 
series  were  then  taken  with  a  time  interval  of  69  seconds.  Each  series  contained  64 
coronal  slices  with  a  matrix  of  128  x256  pixels  and  an  in-plane  resolution  of 
1.25mmxl.25mm.  Slice  thickness  ranged  from  2  to  3  mm  depending  on  breast  size.  For 
this  database  121  primary  mass  lesions  have  been  outlined  by  an  experienced  radiologist, 
77  of  the  lesions  are  malignant  and  44  lesions  are  benign,  as  revealed  by  biopsy. 

We  are  also  collecting  more  data  in  the  University  of  Chicago  Hospitals.  The  imaging 
technique  is  the  same  as  Database-I  above.  We  have  collected  1 80  patients’  images  in 
the  past  year  and  more  data  are  being  collected  currently. 

2.  Investigation  of  image  features 

We  used  Database-I  for  investigation  of  features  that  contribute  to  the  diagnosis  of 
breast  cancer.  The  lesions  were  delineated  by  an  experienced  radiologist  as  well  as 
independently  by  computer  using  an  automatic  volume-growing  algorithm.  Fourteen 
features  that  were  extracted  automatically  from  the  lesions  could  be  grouped  into  three 
categories  based  on:  (I)  morphology;  (II)  enhancement  kinetics;  and  (III)  time  course  of 
enhancement-variation  over  the  lesion.  A  stepwise  feature  selection  procedure  was 
employed  to  select  an  effective  subset  of  features,  which  were  then  combined  by  linear 
discriminant  analysis  (LDA)  into  a  discriminant  score,  related  to  the  likelihood  of 
malignancy.  The  classification  performances  of  individual  features  and  the  combined 
discriminant  score  were  evaluated  with  receiver  operating  characteristic  (ROC)  analysis. 
With  the  radiologist-delineated  lesion  contours,  stepwise  feature  selection  yielded  4 
features  and  an  Az  value  of  0.80  for  the  LDA  in  leave-one-out  cross-validation  testing. 
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With  the  computer-segmented  lesion  volumes,  it  yielded  6  features  and  an  Az  value  of 
0.86  for  the  LDA  in  the  leave-one-out  testing.  A  full  description  of  this  study  is  given  in 
reference  [5]  which  is  attached  as  Appendix  A. 

3.  Automatic  identification  of  characteristic  enhancement  kinetics  in  breast  MR 
lesions 

In  [5]  we  have  demonstrated  the  usefulness  of  enhancement  kinetics  features  in 
distinguishing  between  malignant  and  benign  lesions.  Due  to  uptake  heterogeneity, 
however,  the  enhancement  curve  obtained  from  a  specific  region  in  the  lesion  may 
outperform  that  from  averaging  over  the  entire  lesion.  We  have  developed  an  automatic 
method,  based  fuzzy  c-means  (FCM)  clustering  analysis,  to  identify  characteristic 
enhancement  kinetics  within  the  lesion.  Our  results  have  shown  that  the  developed 
methods  significantly  improved  the  performance  of  enhancement  kinetics  in  the  diagnosis 
of  breast  cancer.  A  more  detailed  summary  of  these  methods  can  be  found  in  reference 
[6]  which  is  also  attached  as  Appendix  B. 

4.  Development  of  methods  for  shading  artifacts  correction  in  breast  MR  images 

Most  clinical  MR  images  are  corrupted  by  shading  artifacts,  which  is  slow  intensity 
variations  of  the  same  tissue  over  the  image  domain.  Such  artifacts  might  affect 
computerized  analysis  such  as  tissue  segmentation.  We  developed  a  FCM  based 
approach  that  simultaneously  estimates  the  shading  effect  while  segmenting  the  image. 
A  full  description  of  the  methods  is  in  reference  [7]  which  is  attached  as  Appendix  C. 
Further  study  will  be  done  to  investigate  the  influence  of  the  shading  artifacts  on  lesion 
extent  assessment  and  on  the  performance  of  features  in  the  task  of  distinguishing 
between  malignant  and  benign  lesions. 

5.  Development  of  computerized  methods  for  assessment  of  tumor  extent 

The  assessment  of  tumor  extent  is  a  key  step  in  computerized  analysis  of  breast 
lesions,  as  the  features  are  derived  from  the  segmented  lesion.  In  our  previous  study  [5], 
we  used  a  region-growing  method  for  the  segmentation  of  breast  lesions.  We  developed 
a  FCM  based  algorithm,  which  assesses  the  tumor  extent  based  on  analysis  of  the  signal 
intensity-curves  within  a  region  of  interest.  The  method  has  been  shown  to  be  more 
accurate  than  the  region-growing  method  in  the  assessment  of  lesion  extent  using 
radiologist  outlined  lesion  as  a  ground  truth.  The  preliminary  study  has  been  presented  at 
the  18th  International  Congress  and  Exhibition  of  Computer  Assisted  Radiology  and 
Surgery,  (CARS  2004),  and  the  paper  [8]  has  been  published  in  the  Proceedings  of  CARS 
2004.  This  short  paper  is  attached  in  Appendix  D. 

A  complete  study  for  this  topic  has  been  submitted  to  journal  IEEE  Transactions  on 
Medical  Imaging  and  the  manuscript  is  under  revision  [9]. 

6.  Development  of  a  prototype  workstation  for  computer-aided  diagnosis  (CAD)  in 
contrast-enhanced  MRI 

We  are  developing  a  prototype  workstation  for  CAD  in  CE-MRI,  which  integrates 
our  algorithms  for  analysis  of  breast  MR  images  with  a  user-friendly  graphical  user 
interface.  The  current  version  of  this  software  package  has  included  image  visualization, 
automatic  lesion  segmentation,  and  automatic  characteristic  kinetic  curve  generation.  It 
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has  been  presented  in  the  InfoRad  Exhibit  in  the  99th  RSNA  conference  and  has  received 
positive  comments  from  both  researchers  and  radiologists. 
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KEY  RESEARCH  ACCOMPLISHMENTS 

•  Collection  and  maintenance  of  clinical  database:  dynamic  MR  images  of  over  200 
patients  have  been  collected  which  are  well  suitable  for  the  proposed  research; 

•  Investigation  of  lesion  features:  a  new  family  of  features  have  been  investigated 
and  proved  to  be  useful  in  the  diagnosis  of  breast  cancer; 

•  Developed  automatic  methods  for  identification  of  characteristic  enhancement 
kinetics  in  breast  MR  lesions  which  outperformed  the  performance  of 
enhancement  kinetics  obtained  from  averaging  over  the  entire  lesion  in  the 
diagnosis  of  breast  cancer; 

•  Developed  automatic  methods  for  shading  artifacts  correction  in  breast  MR 
images  which  will  potentially  increase  the  accuracy  of  lesion  segmentation  and 
the  performance  of  lesion  features  for  diagnosis; 

•  Developed  computerized  methods  for  assessment  of  tumor  extent  which 
outperformed  the  previous  region-growing  method. 
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REPORTABLE  OUTCOMES 

Peer  reviewed  journal  papers 

•  W.  Chen,  M.  L.  Giger,  L.  Lan,  and  U.  Bick,  “Computerized  interpretation  of 
breast  MRI:  Investigation  of  enhancement-variance  dynamics,”  Medical  Physics , 
3 1(5):  1076- 1082,  (2004). 

•  W.  Chen,  M.  L.  Giger,  U.  Bick,  “A  fuzzy  c-means  (FCM)  based  approach  for 
computerized  segmentation  of  breast  lesions  in  contrast-enhanced  MR  images”, 
IEEE  Trans.  On  Med.  Imag.,  in  revision. 

Papers  in  conference  proceedings 

•  W.  Chen,  M.  L.  Giger,  and  U.  Bick,  “Automated  identification  of  temporal 
pattern  with  high  initial  enhancement  in  dynamic  MR  breast  lesions  using  fuzzy 
c-means  algorithm,”  Proceedings  SPIE,  2004. 

•  W.  Chen,  M.  L.  Giger,  “A  fuzzy  c-means  (FCM)  based  algorithm  for  intensity 
inhomogeneity  correction  and  segmentation  of  MR  images,”  Proceedings  of  IEEE 
International  Symposium  on  Biomedical  Imaging,  2004. 

•  W.Chen,  M.L. Giger,  G.  Newstead,  U.  Bick,  L.  Lan,  “Computerized  assessment  of 
tumor  extent  in  contrast-enhanced  MR  images  of  the  breast,”  Proceedings  of  18th 
International  Congress  and  Exhibition-  Computer  Assisted  Radiology  and 
Surgery,  (CARS  2004),  in  press 

Presentations 

•  W.  Chen,  M.  L.  Giger,  and  U.  Bick,  “Automated  identification  of  temporal 
pattern  with  high  initial  enhancement  in  dynamic  MR  breast  lesions  using  fuzzy 
c-means  algorithm,”  SPIE-Medical  Imaging,  San  Diego,  California,  February 
2004. 

•  W.  Chen,  M.  L.  Giger,  G.  M.  Newstead,  U.  Bick,  L.  Lan,  “Computerized  analysis 
of  contrast-enhanced  MR  images  of  the  breast:  automated  bias  field  correction 
and  identification  of  characteristic  signal-time  curves,”  90th  Assembly  and 
Annual  Meeting  of  Radiological  Society  of  North  America,  Chicago,  Illinois, 
USA,  December  2004. 

Posters/Education  Exhibits 

•  W.  Chen,  M.  L.  Giger,  “A  fuzzy  c-means  (FCM)  based  algorithm  for  intensity 
inhomogeneity  correction  and  segmentation  of  MR  images,”  IEEE  International 
Symposium  on  Biomedical  Imaging,  Arlington,  Virginia,  April  2004. 

•  W.Chen,  M.L. Giger,  G.  Newstead,  U.  Bick,  L.  Lan,  “Computerized  assessment  of 
tumor  extent  in  contrast-enhanced  MR  images  of  the  breast,”  18th  International 
Congress  and  Exhibition-  Computer  Assisted  Radiology  and  Surgery  (CARS 
2004),  Chicago,  Illinois,  June  2004. 
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•  W.Chen,  M.L.Giger,  G.  Newstead,  L.  Lan  ,  “Computer-aided  aiagnosis  of  breast 
cancer  with  dynamic  contrast-enhanced  MRI”,  90™  Assembly  and  Annual 
Meeting  of  Radiological  Society  of  North  America,  Chicago,  Illinois,  USA, 
December  2004. 


Awards 


•  Doolittle-Harrison  Fellowship  ($500),  The  University  of  Chicago 
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CONCLUSIONS 

The  recipient  of  the  Predoctoral  Traineeship  Award  has  taken  all  the  required  core 
courses  and  many  research  related  elective  courses  as  well.  These  trainings  have  proven 
useful  for  the  recipient  to  achieve  the  proposed  research  goals.  The  goal  of  this  grant  is  to 
support  the  recipient  for  research  towards  Ph.D.  degree.  The  trainee  has  made  good 
progress  in  this  aspect,  as  he  has  written  up  his  dissertation  proposal,  which  is  a  key  step 
to  obtain  Ph.D.  degree  in  medical  physics  at  the  University  of  Chicago. 

During  the  first  two  years,  we  have  collected  and  maintained  a  dynamic  breast  MRI 
database  that  is  well  suitable  for  the  proposed  research  on  computerized  interpretation  of 
breast  MR  images.  We  have  developed  computerized  methods  for  correction  of  shading 
artifacts,  tumor  segmentation,  and  feature  extraction.  The  results  have  shown  that 
computerized  analysis  and  interpretation  methods  have  great  promise  in  increasing  the 
objectivity,  efficiency,  and  accuracy  of  the  diagnosis  of  breast  cancer. 

Overall,  we  have  well  achieved  the  goals  of  the  first  two  years  and  laid  a  strong 
foundation  for  the  research  in  the  next  years.  Our  goals  in  the  next  year  include 
collection  of  more  data,  investigation  of  features  that  would  better  distinguish  between 
benign  and  malignant  lesions,  and  evaluation  of  the  overall  systematic  methods  with  an 
even  larger  database. 
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Computerized  interpretation  of  breast  MRI:  investigation  of  enhancement- 
variance  dynamicsa) 

Weijie  Chen,  Maryellen  L.  Giger,b)  Li  Lan,  and  Ulrich  Bick 

Department  of  Radiology,  The  University  of  Chicago,  MC  2026,  5841  South  Maryland  Avenue,  Chicago, 

Illinois  60637 

(Received  23  July  2003;  revised  12  February  2004;  accepted  for  publication  13  February  2004; 
published  8  April  2004) 

The  advantages  of  breast  MRI  using  contrast  agent  Gd-DTPAin  the  diagnosis  of  breast  cancer  have 
-  been  well  established.  The  variation  of  interpretation  criteria  and  absence  of  interpretation  guide¬ 
lines,  however,  is  a  major  obstacle  for  applications  of  MRI  in  the  routine  clinical  practice  of  breast 
imaging.  Our  study  aims  to  increase  the  objectivity  and  reproducibility  of  breast  MRI  interpretation 
by  developing  an  automated  interpretation  approach  for  ultimate  use  in  computer-aided  diagnosis. 

The  database  in  this  study  contains  121  cases:  77  malignant  and  44  benign  masses  as  revealed  by 
biopsy.  Images  were  obtained  using  a  Tl-weighted  3D  spoiled  gradient  echo  sequence.  After  the 
acquisition  of  the  precontrast  series,  Gd-DTPA  contrast  agent  was  injected  intravenously  by  power 
injection  with  a  dose  of  0.2  mmol/kg.  Five  postcontrast  series  were  then  taken  with  a  time  interval 
of  60  s.  Each  series  contained  64  coronal  slices  with  a  matrix  of  128X256  pixels  and  an  in-plane 
resolution  of  1. 25 X  1.25  mm2.  Slice  thickness  ranged  from  2  to  3  mm  depending  on  breast  size. 

The  lesions  were  delineated  by  an  experienced  radiologist  as  well  as  independently  by  computer 
using  an  automatic  volume-growing  algorithm.  Fourteen  features  that  were  extracted  automatically 
from  the  lesions  could  be  grouped  into  three  categories  based  on  (I)  morphology,  (II)  enhancement 
kinetics,  and  (HI)  time  course  of  enhancement-variation  over  the  lesion.  A  stepwise  feature  selec¬ 
tion  procedure  was  employed  to  select  an  effective  subset  of  features,  which  were  then  combined  by 
linear  discriminant  analysis  (LDA)  into  a  discriminant  score,  related  to  the  likelihood  of  malig¬ 
nancy.  The  classification  performances  of  individual  features  and  the  combined  discriminant  score 
were  evaluated  with  receiver  operating  characteristic  (ROC)  analysis.  With  the  radiologist- 
delineated  lesion  contours,  stepwise  feature  selection  yielded  four  features  and  an  A.  value  of  0.80 
for  the  LDA  in  leave-one-out  cross-validation  testing.  With  the  computer-segmented  lesion  vol¬ 
umes,  it  yielded  six  features  and  an  A-  value  of  0.86  for  the  LDA  in  the  leave-one-out 
testing.  ©  2004  American  Association  of  Physicists  in  Medicine.  [DOI:  10.1118/1.1695652] 

Key  words:  breast  cancer,  contrast-enhanced  MRI,  Gd-DTPA,  computer-aided  diagnosis  (CAD), 

ROC  analysis 


I.  INTRODUCTION 

Breast  MRI  has  emerged  as  a  promising  modality  for  the 
detection  and  diagnosis  of  breast  cancer  since  the  introduc¬ 
tion  of  gadolinium-diethylenetriamine  penta-acetic  acid  (Gd- 
DTPA)  as  a  contrast  agent. 1-3  Contrast-enhanced  MRI  (CE- 
MRI)  allows  lesions  to  be  distinguished  from  normal  tissues 
due  to  the  increased  vascularity  and  capillary  permeability  of 
tumors.  CE-MRI  offers  three-dimensional  spatial  informa¬ 
tion  and  temporal  information  of  breast  cancer,  qualifying  it 
as  an  encouraging  complementary  modality  to  conventional 
imaging  methods,  such  as  x-ray  mammography  and  sonog¬ 
raphy. 

Despite  its  well-recognized  advantages,  applications  of 
MRI  in  the  routine  clinical  practice  of  breast  imaging  are 
limited.  One  of  the  most  important  obstacles  is  the  lack  of 
interpretation  guidelines;4,5  very  few  attempts  have  been 
made  to  standardize  the  interpretation  of  breast  MR  images. 
Among  the  few  efforts  in  this  regard  is  the  work  of  the  In¬ 
ternational  Working  Group  on  Breast  MRI,  which  has  devel¬ 
oped  and  validated  a  detailed  lexicon  for  breast  MRI 
interpretation.6  Investigators  use  a  large  variety  of  diagnostic 


criteria2, 3,7-13  that  help  classify  lesions  as  benign  or  malig¬ 
nant.  The  interpretation  criteria  in  the  current  literature  fall 
into  two  major  categories:5,14  morphologic  features8,9,13  and 
enhancement  kinetics,2,3,12  i.e.,  the  time  course  of  signal  in¬ 
tensity  within  the  suspected  lesions.  While  the  studies  using 
these  criteria  have  shown  promising  results  (sensitivity  from 
92%  to  100%,  specificity  from  53%  to  100%)  as  reported  in 
a  recent  review,15  significant  variation  does  exist.  The  varia¬ 
tion  of  the  results  may  be  due  to  two  reasons  from  the  inter¬ 
pretation  aspect;  one  is  the  interobserver  variation  from  dif¬ 
ferent  subjective  judgments,16-18  and  the  other  is  that  the 
current  interpretation  schemes  might  not  be  sufficiently  ro¬ 
bust. 

The  aim  of  computerized  interpretation  of  medical  images 
is  to  obtain  quantitative  indices  of  malignancy.  It  has  the 
advantage  of  being  objective,  automatic,  and,  furthermore,  it 
may  provide  unique  information  that  might  be  difficult  to 
assess  visually,  especially  for  time-series  3D  MR  images. 

This  study  aims  to  use  computerized  methods  to  investi¬ 
gate  the  potential  of  enhancement-variance  dynamics  in  the 
interpretation  of  contrast-enhanced  breast  MR  images.  Auto- 
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matically  extracted  features  are  based  on  the  time  course  of 
enhancement-variance  within  lesions.  Features  based  on  en¬ 
hancement  kinetics  and  morphology  of  lesions  are  studied  as 
well  for  comparison.  Also  compared  are  the  performances  of 
features  extracted  from  lesions  delineated  by  experienced  ra¬ 
diologist  and  those  from  lesions  segmented  by  the  computer 
using  a  volume-growing  algorithm.  Finally,  different  features 
are  selected  and  merged  into  an  estimate  of  malignancy  us¬ 
ing  automated  classification. 

II.  MATERIALS  AND  METHODS 

A.  Image  database 

The  database  in  this  study  contains  121  cases:  77  malig¬ 
nant  and  44  benign  masses  as  revealed  by  biopsy.  Images 
were  obtained  using  a  T1 -weighted  3D  spoiled  gradient  echo 
sequence  (TR=8.1ms,  TE=4ms,  flip  angle=30°).  Fat 
suppression  was  not  employed.  The  patients  were  scanned  in 
prone  position  using  a  standard  double  breast  coil  on  a  1 .5  T 
whole-body  MRI  system  (Siemens  Vision,  Siemens,  Erlan¬ 
gen,  Germany).  After  the  acquisition  of  the  precontrast  se¬ 
ries,  Gd-DTPA  contrast  agent  was  delivered  intravenously  by 
power  injection  with  a  dose  of  0.2  mmol/kg  and  a  flow  rate 
of  2  ml/s.  Injection  of  contrast  was  followed  by  a  saline  flush 
of  20  ml  with  the  same  flow  rate.  Five  postcontrast  series 
were  then  taken  with  a  time  interval  of  60  s.  Each  series 
contained  64  coronal  slices  with  a  matrix  of  128 
X  256  pixels  and  an  in-plane  resolution  of  1 .25  X  1 .25  mm2 . 
Slice  thickness  ranged  from  2  to  3  mm  depending  on  breast 
size.  The  image  database  had  been  retrospectively  collected 
under  an  IRB-approved  protocol. 

B.  Methods 

The  computerized  interpretation  scheme  used  in  this  study 
begins  with  the  segmentation  of  the  lesion  within  the  image. 
The  suspect  masses  were  delineated  both  manually  by  an 
experienced  radiologist  (U.B.)  and  automatically  by  the  com¬ 
puter  using  a  3D  volume-growing  algorithm.19  Next,  mul¬ 
tiple  features  that  characterize  the  spatial  and  kinetic  proper¬ 
ties  of  the  lesions,  and  thus  would  potentially  help 
differentiate  the  malignant  cases  from  the  benign  cases,  were 
extracted  automatically.  Stepwise  feature  selection20  was  em¬ 
ployed  to  select  a  set  of  features  that  perform  efficiently  in 
classifying  the  lesions  as  malignant  or  benign.  Using  linear 
discriminant  analysis21  (LDA),  the  selected  features  were 
then  merged  into  a  single  numerical  value  that  is  related  to 
estimated  likelihood  of  malignancy. 

We  evaluated  the  classification  performance  of  individual 
features  and  the  merged  discriminant  score  in  the  task  of 
distinguishing  between  malignant  and  benign  lesions  by  us¬ 
ing  receiver  operating  characteristics  (ROC)  analysis.22  The 
areaunder  the  maximum  likelihood-fitted  binormal  ROC 
curve,  A, ,  was  used  as  an  index  of  performance.23  We  used 
the  CLABROC24,25  algorithm  to  determine  the  statistical  sig¬ 
nificance  (p  value)  of  the  difference  between  two  A.  values, 
associated  with  two  lesion  delineation  methods.  The  CLA- 


Table  I.  Features  investigated  in  this  study  are  grouped  into  three  major 
categories. 


I.  Morphologic 

Fj  i :  Maximum  std.  of  RGH  value 
Circularity 
F,  j :  Irregularity 
F,y.  Margin  gradient 
Fi  5  :  Variance  of  margin  gradient 

II.  Enhancement  kinetics 
F„,  i :  Maximum  uptake 
F,i2-  Peak  location 
F„,  j :  Uptake  rate 

F II, 4  ■  Washout  rate 

III.  Enhancement-variance  dynamics 

Pm ,i  :  Maximum  enhancement-variance 

Fin,  2 :  Peak  location 

F i„.i  :  Increasing  rate 

Fin,  4  •'  Decreasing  rate 

Fin,  5  :  Enhancement- variance  at  time  #1 


BROC  algorithm  uses  a  univariate  z-score  test  to  test  the 
difference  between  areas  under  two  ROC  curves. 

1.  Lesion  delineation 

For  the  manual  delineation,  a  radiologist  (U.B.),  blinded 
to  the  histological  diagnosis,  contoured  the  enhanced  tumor 
area  in  each  slice  that  intersected  the  lesion  in  the  subtraction 
images  (postcontrast  minus  precontrast),  using  the  nonsub- 
tracted  MR  images  as  reference. 

For  the  automatic  segmentation  by  computer,  we  used  a 
volume-growing  based  algorithm.19  In  this  algorithm,  the 
breast  volume  is  first  segmented  at  a  threshold  derived  from 
the  global  histogram  of  voxel-values  by  maximizing  the  in¬ 
terclass  variance  between  breast  and  background.  Then  the 
border  of  the  segmented  breast  volume  is  removed  by  mor¬ 
phological  erosion  using  a  3X3X3  structuring  element. 
Next,  regions  with  high  contrast  uptake  are  enhanced  by 
voxelwise  computation  of  the  variance  of  the  voxel  values 
over  time.  A  spherical  region  of  interest  (SROI)  that  encom¬ 
passes  the  enhanced  region  is  then  automatically  selected 
from  a  set  of  spheres  expanding  outwards  from  a  manually 
selected  seed  point.  From  the  preprocessed  volume  within 
the  SROI,  a  segmentation  threshold  is  then  computed  by 
maximizing  the  interclass  variance  between  enhanced  lesion 
and  background  voxels  in  the  SROI.  Finally,  6-point- 
connected  volume  growing  is  performed  to  yield  the  3D  seg¬ 
mented  lesion. 

2.  Computerized  feature  extraction 

In  this  study,  our  primary  interest  is  to  investigate  the 
potential  of  enhancement  variance  dynamics  features  in  the 
classification  of  suspicious  lesions  as  malignant  or  benign. 
Other  features  are  studied  as  well  for  comparison.  Features 
investigated  in  this  study,  as  listed  in  Table  I,  are  grouped 
into  three  categories:  (I)  morphological  features,  (II)  en¬ 
hancement  kinetics  based  features,  and  (IE)  features  related 
to  the  time  course  of  enhancement-variance  over  the  lesion. 
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Fig.  1 .  Examples  of  lesion  segmenta¬ 
tion  and  the  corresponding  enhance¬ 
ment  kinetics  and  enhancement- 
variance  dynamics  curves:  (a) 
malignant  case  and  (b)  benign  case.  In 
each  case,  a  central  slice  image  of  the 
volume  lesion  is  shown  (upper  left). 
Also  shown  are  the  radiologist- 
outlined  lesions  and  the  corresponding 
lesion  enhancement  [Eq.  (1)]  at  six 
time  points  (middle  left),  the 
computer-segmented  lesion,  and  the 
corresponding  lesion  enhancement  at 
six  time  points  (lower  left).  The  asso¬ 
ciated  kinetics  curves  are  shown 
(right). 


The  features  in  the  first  category  have  been  described  in  de¬ 
tail  in  a  previous  paper  from  our  group12  and  are  only  sum¬ 
marized  here. 

(i)  Morphological  features .n  The  feature  “maximum 
standard  deviation  (std)  of  radial  gradient  histogram 
(RGH)”  (Fi  x)  quantifies  how  well  the  image  structures  in  a 
suspected  lesion  extend  in  a  spherical  pattern  originating 
from  the  center  of  the  lesion.  The  feature  “circularity” 
(Fj  2)  measures  conformity  of  a  lesion  to  a  spherical  shape 
and  the  feature  “irregularity”  (F,f)  indicates  the  roughness 
of  the  lesion  surface.  The  features  “margin  gradient”  (FI  4) 
and  “variance  of  margin  gradient”  (FI  5)  are  related  to  the 
sharpness  of  the  lesion  margin  by  evaluating  voxel-value 
gradients  and  their  variations  along  the  margin  of  the  suspect 
lesion. 

(ii)  Enhancement  kinetics  based  features.  Enhancement 
kinetics  is  related  to  the  time  course  of  signal  intensity 
within  the  lesion.  Denote  S(r,i)  as  the  voxel  value  at  loca¬ 
tion  r  in  the  lesion  at  time  frame  i,  i  runs  from  0  (i.e.,  the 
precontrast  frame)  to  5  (i.e.,  the  last  postcontrast  frame).  For 


each  voxel  in  the  lesion,  the  contrast  enhancement  is  com¬ 
puted: 


C(r,z)  = 


S(r,z)-S(r,0) 
S(  r,0)  ’ 


/  =  0,1,...,  5. 


(1) 


C(r,z)  is  a  quantity  that  has  been  shown26  to  be  related  to 
Gd-DTPA  concentration  in  the  extracellular  space  of  breast 
tissue  at  voxel  r.  In  special  circumstances  this  relation  ap¬ 
proaches  linearity.  Note  that  at  i  =  0,C(r,z)  =  0.  The  enhance¬ 
ment  dynamics  can  be  described  by  the  average  enhance¬ 
ment  over  the  lesion  at  each  time  point,  i.e., 

1  1 

C(/)=7S  C(r,z),  z  =  0,l,. ..,5,  (2) 

L  r=  1 

where  L  is  the  number  of  voxels  in  the  segmented  lesion. 

Four  features  are  derived  from  the  enhancement  kinetics. 
The  maximum  uptake  (Fn  l)  is  the  maximum  enhance- 
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ment,  i.e.,  Fnl  =  max/=015C(/).  The  time  frame  index  at 
which  the  maximum  enhancement  occurs  is  a  feature  and  is 
referred  to  here  as  “peak  location”  (F//2). 

The  uptake  rate  (Fni)  of  the  contrast  agent  is  defined  as 

Fu,3=Fju/F„'2.  (3) 

The  washout  rate  ( FlI  4 )  of  the  contrast  agent  is  defined 
as 


77,4 ' 


Fty-CjS) 

5  ”P//,2 

0 


(P//,2*5), 


(^77,2=5). 


(4) 


(in)  Enhancement-variance  dynamics  features. 
Enhancement-variance  dynamics  describes  the  time  course 
of  the  spatial  variance  of  the  enhancement  within  the  lesion 
and  is  defined  by 

1  L 

^)=rr2  [C(r,i)-c(/)]2,  ;=o,i,...,5.  (5) 

E  —  1  r=l 


Five  features  are  derived  from  the  enhancement-variance 
dynamics.  The  maximum  variation  of  enhancement 
(F/// ,)  is  the  maximum  spatial  variance  of  enhancement, 
i.e.,  Fm  \  =  max,- =0,1.,., 5^(0-  The  time  frame  index  at  which 
the  maximum  variance  occurs  is  a  feature  and  is  referred  to 
here  as  “peak  location”  (Fm2)  of  the  enhancement- 
variance  dynamics. 

The  enhancement-variance  increasing  rate  (Fmf)  de¬ 
scribes  how  fast  the  enhancement-variation  within  the  lesion 
reaches  maximum,  defined  by 

Fm, i = Fin ,  i  /  F ui2 .  (6) 
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Fio.  2.  Performance  of  the  computer  segmentation  as  compared  to  radiolo¬ 
gist  segmentation.  Curve  indicates  the  fraction  of  lesions  correctly  seg¬ 
mented  as  a  function  of  the  overlap  criterion. 


The  enhancement-variance  decreasing  rate  (F///4)  in¬ 
dicates  how  fast  the  enhancement-variance  decreases  from 
the  maximum,  defined  by 


■P/7/,4” 


P///,i”P(5) 
^  F  „1,2 


(Fm, 2^5), 


(Fi„,  2=5). 


(7) 


Table  II.  Performance  of  14  computer-extracted  features  in  distinguishing  between  malignant  and  benign  breast 
lesions  that  are  delineated  by  both  radiologist  (column  2)  and  computer  (column  3).  The  value  after  “±”  is  the 
standard  deviation  (s.d.)  associated  with  each  A-  value.  The  two  tailed  /> -value  was  calculated  using  the 
univariate  z -score  test. 


Feature 

A.±  1  s.d.: 
Radiologist 
outlined 

A-±  I  s.d: 
Computer 
segmented 

p -value 

I.  Morphologic 

Fr] :  Maximum  std.  of  RGH  value 

0.59+0.05 

0.56±0.05 

0.33 

F/  2 :  Circularity 

0.57±0.06 

0.65±0.05 

0.42 

f/j :  Irregularity 

0.66±0.05 

0.54±0.06 

0.19 

f'tA :  Margin  gradient 

0.58±0.06 

0.60±0.06 

0.57 

Ft, 5 :  Variance  of  margin  gradient 

0.60+0.05 

0.51  ±0.06 

0.03 

II.  Enhancement  kinetics 

Fn  i :  Maximum  uptake 

0.55±0.06 

0.51  ±0.06 

0.25 

F n  2  '■  Peak  location 

0.75+0.05 

0.79±0.05 

0.15 

F//j :  Uptake  rate 

0.68  ±0.05 

0.66±0.05 

0.50 

F,i  4 :  Washout  rate 

0.73  ±0.05 

0.78±0.05 

0.09 

III.  Enhancement-variance  dynamics 

Fm  i :  Maximum  variation  of  enhancement 

0.50±0.06 

0.52  ±0.06 

0.86 

Fm, 2 :  Peak  location 

0.74±0.05 

0.77  ±0.05 

0.08 

F,n  j :  Increasing  rate 

0.58±0.06 

0.58  ±0.06 

0.96 

Fm  4 :  Decreasing  rate 

0.74±0.05 

0.73  ±0.05 

0.90 

F m i5 :  Enhancement-variance  at  time  #1 

0.65  ±0.06 

0.63  ±0.06 

0.31 
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Fio.  3.  Performance  of  the  14  computer-extracted  features  in  distinguishing 
between  malignant  and  benign  breast  lesions  that  are  delineated  by  both 
radiologist  (solid)  and  computer  (dash). 


The  enhancement-variance  at  the  first  postcontrast 
frame  (/r//7 5),  V(l),  reveals  the  uptake  inhomogeneity  at 
the  early  phase  of  uptake. 

Figure  1  shows  two  examples  (one  malignant  and  one 
benign)  of  lesion  segmentation  and  the  corresponding  en¬ 
hancement  kinetics  curves  and  enhancement-variance  curves. 
Note  that  the  lesion  sizes  in  the  two  examples  are  different, 
so  in  the  figure  the  pixels  of  the  smaller  lesion  [Fig.  1(b)] 
look  larger. 

III.  RESULTS 

Performance  of  the  computer  segmentation  is  shown  in 
Fig.  2.  Here  overlap  is  defined  as  the  volume  of  intersection 
of  the  radiologist-delineated  lesion  and  the  computer- 
segmented  lesion  divided  by  the  volume  of  their  union.  In 
this  database  with  121  mass  lesions,  84%  of  the  lesions  were 
correctly  segmented  at  an  overlap  threshold  of  0.4.  The  two 
segmentation  examples  shown  in  Fig.  1  had  overlaps  of  (a) 
0.77  and  (b)  0.17,  respectively. 

Table  II  and  Fig.  3  show  the  A.  values  and  the  associated 
standard  deviations  indicating  individual  performance  levels 
of  the  14  features  in  the  task  of  distinguishing  between  ma¬ 
lignant  and  benign  breast  lesions  that  were  delineated  by 
either  radiologist  or  computer.  The  results  indicate  that  all 
three  categories  of  features  show  potential  for  the  classifica¬ 
tion  task. 

As  also  demonstrated  in  Table  II  and  Fig.  3,  at  the  signifi¬ 
cance  level  p  —  0.05,  we  failed  to  show  a  statistically  signifi- 
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cant  difference  between  the  performance  (i.e.,  A.)  of  features 
obtained  using  the  radiologist-outlined  lesions  and  the  per¬ 
formance  of  those  obtained  with  the  computer-segmentation 
(except  for  feature  F,  5)  in  the  task  of  distinguishing  between 
malignant  and  benign  breast  lesions. 

Stepwise  feature  selection20  selected  two  sets  of 
features — one  set  for  each  of  the  two  methods  of  lesion  de¬ 
lineation  (Table  HI).  The  selected  feature-set  from  features 
based  on  the  radiologist-outlined  lesions  (set  1)  includes  four 
features:  peak  location  of  enhancement-variance  dynamics 
(Fmi),  irregularity  (F/  3),  washout  rate  of  enhancement  ki¬ 
netics  (Fn  4),  and  peak  location  of  enhancement  kinetics 
(Fn 2).  The  leave-one-out  cross-validation  using  linear  dis¬ 
criminant  analysis  to  merge  the  selected  features  yields  an  A, 
value  of  0.80  in  the  task  of  distinguishing  between  malignant 
and  benign  lesions.  The  selected  feature-set  from  features 
based  on  the  computer-segmented  lesions  (set  2)  includes  six 
features:  peak  location  of  enhancement-variance  dynamics 
(F  111,2),  enhancement-variance  increasing  rate  (F!ni),  peak 
location  of  enhancement  kinetics  ( F„2 ),  margin  gradient 
(Fj 4),  variance  of  margin  gradient  (FI5),  and  maximum 
standard  deviation  of  RGH  value  (F/ ,).  The  leave-one-out 
cross-validation  using  linear  discriminant  analysis  to  merge 
the  selected  features  yields  an  A.  value  of  0.86. 

IV.  DISCUSSION 

Enhancement  kinetic  analysis  evaluates  how  the  contrast 
enhancement  within  the  lesion  changes  in  a  period  of  time, 
and  reveals  the  uptake  and  washout  characteristics  of  the 
intravenous  contrast  within  a  lesion.  Studies  have  shown  that 
malignant  cases  tend  to  have  rapid  uptake  and  washout.12 
Our  results  are  consistent  with  such  published  results.  In  our 
analysis,  64%  of  the  benign  lesions  showed  no  washout  over 
the  six  acquisition  frames,  i.e.,  their  enhancement-curves 
keep  increasing.  For  malignant  lesions  only  24%  demon¬ 
strated  no  washout.  It  is  worth  noting  that  all  the  features  are 
calculated  based  on  the  average  enhancement  of  the  entire 
lesion.  In  clinical  practice,  however,  radiologists  can  choose 
any  region  with  which  to  assess  the  enhancement  kinetics. 
Our  future  work  of  interest  is  to  develop  automatic  methods 
to  extract  region  of  maximum  enhancement.  The 
enhancement-kinetics  based  features  on  such  “hot  areas”  are 
expected  to  perform  better  in  differential  diagnosis. 

Enhancement-variance  dynamics  reveals  how  the  spatial 
distribution  of  contrast  enhancement  in  the  lesion  region 
changes  in  a  period  of  time.  Our  results  showed  that  the 
spatial  variance  of  enhancement  peaks  earlier  in  malignant 
lesions  (Fig.  4).  The  enhancement-variance  in  a  particular 


Table  III.  Stepwise  feature  selection  and  the  performance  of  the  merged  features  using  leave-one-out  cross- 
validation.  The  value  after  “±”  is  the  standard  deviation  (s.d.)  associated  with  each  A.  value.  The  two-tailed 


p -value  was  calculated  using  the  univariate  2 -score  test. 

Radiologist  outlined  lesions 

Computer  segmented  lesions 

Feature  Set  1  A.±  1  s.d. 

^111,2^1,^11,4^11,2  0.80±0.04 

Feature  Set  2  A.±  1  s.d. 

?  111,2? in,y? n.2?  1,4?  1,5?  1, 1  0.86±  0.04 

p -value 
0.07 
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Fig.  4.  Histogram  of  peak  location  of  enhancement-variance  dynamics  for 
the  computer-segmented  lesions. 


time  frame  reveals  the  degree  of  inhomogeneity  of  the  en¬ 
hancement.  Our  results  show  that  malignant  lesions  tend  to 
be  less  homogenous  at  the  first  postcontrast  time  frame 
( Fm,s )•  However,  malignant  and  benign  lesions  showed 
similar  results  for  the  maximum  inhomogeneity  (F,„  [),  as 
determined  over  all  time  frames. 

The  relative  performance  of  the  morphological  and  ki¬ 
netic  features  depends  on  the  MR  imaging  protocol  involved. 
Basically,  morphologic  features  perform  better  in  high  spatial 
resolution  images  and  enhancement  kinetics  features  perform 
better  in  high  temporal  resolution  images.27  Intuitively,  mor¬ 
phological  features  might  be  more  sensitive  to  lesion  delin¬ 
eation  methods.  With  the  imaging  protocol  and  lesion  seg¬ 
mentation  methods  used  in  this  study,  overall,  the  kinetics 
features  appeared  to  perform  better  than  did  the  morphologi¬ 
cal  features  in  the  task  of  distinguishing  between  malignant 
and  benign  lesions. 

At  present,  protocols  used  in  breast  MR  imaging  world¬ 
wide  are  far  from  reaching  a  consensus  and  optimal  imaging 
parameters  are  still  under  investigation.  The  imaging  tech¬ 
nique  used  in  our  study  is  similar  to  that  launched  by  Kaiser 
etal.,2  which  was  called  “the  archetype  of  dynamic  breast 
MRI”  by  Kuhl  et  al.  :5  acquisition  of  one  precontrast  and  a 
series  of  postcontrast  images  including  both  breasts  at  a  tem¬ 
poral  resolution  of  60  s.  A  potential  limitation  of  the  tech¬ 
nique  is  that  the  data  were  acquired  over  only  5  min  after 
contrast  injection,  as  there  is  probably  some  useful  informa¬ 
tion  beyond  5  min.28  The  choice  of  the  parameters  is  a  trade¬ 
off  to  allow  for  a  shorter  exam  time,  which  is  better  accepted 
by  patients  and  is  less  costly.  This  limited  acquisition  time 
problem  might  have  some  influence  on  the  performance  of 
washout  related  features  (F//  4  The  influence  on  the 

performance  of  other  temporal  features  should  be  much  less 
severe  as,  in  most  cases,  the  signal  reaches  plateau  within  5 
min.  We  regard  as  an  interesting  research  question  the  effect 
of  acquisition  time  on  the  performance  of  computer- 
extracted  features  in  distinguishing  between  malignant  and 
benign  lesions. 


V.  CONCLUSION 

Our  investigation  with  automated  feature  extraction  and 
classification  indicates  that  spatial-variance  dynamics  is  a 
promising  family  of  features  in  the  task  of  distinguishing 
between  malignant  and  benign  MR  breast  lesions.  Combin¬ 
ing  the  morphological  features  of  contrast  enhancements,  en¬ 
hancement  kinetics  features,  and  the  enhancement-variance 
dynamics  features  using  computerized  methods  has  the  po¬ 
tential  to  complement  the  interpretation  of  radiologists  in  a 
consistent,  objective,  and  accurate  way. 
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ABSTRACT 

In  contrast-enhanced  (CE)  MRI  of  the  breast,  signal-intensity  time  curves  have  been  proven  useful  in  differen¬ 
tiating  between  benign  and  malignant  lesions.  Due  to  uptake  heterogeneity  in  the  breast  lesion,  however,  the 
signal-intensity  time  curve  obtained  from  a  specific  region  in  the  lesion  may  outperform  that  from  the  entire 
lesion.  In  this  study,  we  propose  the  use  of  fuzzy  c-means  (FCM)  clustering  algorithms  to  reveal  different  tem¬ 
poral  patterns  within  the  breast  lesion.  The  algorithm  finds  fuzzy  cluster  centers  (i.e. ,  temporal  patterns)  and 
assigns  membership  values  to  each  voxel.  The  temporal  pattern  with  maximum  initial  enhancement  is  selected 
as  the  representative  curve  of  the  lesion  and  the  thresholded  membership  map  is  the  identified  region  of  fast 
enhancement.  The  approach  was  applied  to  the  analysis  of  121  lesions  (77  malignant  and  44  benign).  The 
resulting  representative  curves  were  classified  with  linear  discriminant  analysis  (LDA).  The  differentiation  per¬ 
formance  of  LDA  output  in  leave-one-out  cross  evaluation  was  assessed  using  receiver  operating  characteristic 
(ROC)  analysis.  Our  results  show  that  the  use  of  FCM  significantly  improved  the  performance  of  signal-intensity 
time  curves  in  the  task  of  distinguishing  between  malignant  and  benign  lesions. 

Keywords:  Fuzzy  c-means,  breast  cancer,  contrast-enhanced  MRI,  tumor  heterogeneity,  computer-aided  diag¬ 
nosis  (CAD) 


1.  INTRODUCTION 

Breast  cancer  is  the  most  common  cancer  and  the  second  leading  cause  of  cancer  death  among  women  in  the 
Western  countries.  Mammography  is  currently  the  primary  imaging  modality  in  diagnostic  and  screening  prac¬ 
tice  and  has  achieved  significant  success,  though  there  is  considerable  limitations:  8  -  25%  of  cancers  are  missed 
and  70  —  80%  of  biopsies  turn  out  to  be  benign.  Contrast-enhanced  magnetic  resonance  imaging  (CE-MRI) 
using  Gd-DTPA  is  a  promising  complementary  modality  in  breast  imaging1,2,3  and  has  become  more  widely 
used  recently  as  a  supplemental  imaging  technique,4  particularly  in  the  radiographically  dense  breast.  CE-MRI 
offers  three-dimensional  spatial  information  and  temporal  information,  and  has  demonstrated  extremely  high 
sensitivity  for  breast  cancer.  The  specificity  of  CE-MRI,  however,  is  low  and  this  fact  has  led  intense  efforts  to 
investigate  distinguishing  characteristics  between  malignant  and  benign  lesions.5,6,7,8,9,10 

In  CE-MRI  of  the  breast,  signal-intensity  time  curves  have  been  proven  useful  in  differentiating  between 
benign  and  malignant  lesions.8,9  Due  to  uptake  heterogeneity  in  the  breast  lesion,  however,  the  signal-intensity 
time  curve  obtained  from  a  specific  region  in  the  lesion  may  outperform  that  from  the  entire  lesion.  In  this  study, 
we  propose  the  use  of  fuzzy  c-means  (FCM)  clusterings  algorithm  to  reveal  different  temporal  patterns  within 
the  breast  lesion. 

*Send  correspondence  to  M.L.G.,  E-mail:  m-giger@uchicago.edu;  Phone:  +1  (773)  834-5099;  Fax:  +1  (773)  702-0371 
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2.  MATERIALS  AND  METHODS 

2.1.  Database 

Our  database  in  this  study  contains  121  cases:  77  malignant  and  44  benign  masses  as  determined  by  biopsy. 
Images  were  obtained  using  a  Tl-weighted  3D  spoiled  gradient  echo  sequence  ( TR  —  8.1ms,  TE  —  4  ms, 
flip  angle  =  30°).  Fat  suppression  was  not  employed.  The  patients  were  scanned  in  prone  position  using  a 
standard  double  breast  coil  on  a  1.5  T  whole-body  MRI  system  (Siemens  Vision,  Siemens,  Erlangen,  Germany). 
After  the  acquisition  of  the  precontrast  series,  Gd-DTPA  contrast  agent  was  delivered  intravenously  by  power 
injection  with  a  dose  of  0.2  mmol /kg  and  a  flow  rate  of  2  m// s.  Injection  of  contrast  was  followed  by  a  saline  flush 
of  20  ml  with  the  same  flow  rate.  Five  postcontrast  series  were  then  taken  with  a  time  interval  of  60  seconds.  Each 
series  contained  64  coronal  slices  with  a  matrix  of  128  x  256  pixels  and  an  in-plane  resolution  of  1.25  x  1.25  mm2. 
Slice  thickness  ranged  from  2  to  3  mm  depending  on  breast  size.  The  image  database  had  been  retrospectively 
collected  under  an  IRB-approved  protocol. 

2.2.  Tumor  segmentation 

The  mass  lesions  were  manually  delineated  in  each  slice  by  an  experienced  radiologist  (U.B.).  The  manual 
segmentation  was  performed  in  the  subtraction  images  (postcontrast— precontrast)  in  which  the  tumor  areas 
were  enhanced  in  each  slice.  All  five  subtraction  images  were  used  for  this  purpose.  In  addition,  the  precontrast 
images  were  used  as  reference,  and  complementary  information  from  other  imaging  modalities  were  used  as  well. 

2.3.  Fuzzy  c- means  clustering  analysis 

Our  aim  is  to  identify  different  temporal  patterns  of  all  the  voxels  within  the  lesion  using  FCM  algorithm.  In 
this  context,  the  standard  FCM  algorithm11  is  an  optimization  problem  for  partitioning  a  lesion  of  N  voxels, 
X  =  {x;}^  (xj  =  {xio,Xn,  ...,*15}),  into  c  classes  (temporal  patterns). 

c  N 

arg  min  {J(U,  V;X)  =  -  vk||2}  (1) 

’  k=l i=l 


subject  to: 

C 

Y2  uki  =  1 

k= 1 

0  <  Uki  <  1  Vfc,  i 

where  U  is  the  partition  matrix  whose  element  Uki  is  the  membership  of  the  ith  voxel  for  fcth  class.  V  is  the  cen¬ 
troid  matrix  whose  kth  row  =  {vko,Vki,  ■■■iVks}  is  the  centroid  vector  of  kth  class,  i.e.,  the  temporal  pattern 
prototype  of  kth  class.  The  parameter  p,  called  fuzzy  index,  is  a  weighting  exponent  on  each  fuzzy  membership 
and  determines  the  amount  of  “fuzziness”  of  the  resulting  partition.  The  norm  operator  ||  ■  ||  represents  the 
standard  Euclidean  distance.  The  objective  function  J  is  minimized  when  high  membership  values  are  assigned 
to  the  voxels  whose  intensities  are  close  to  the  centroid  of  its  particular  class,  and  low  membership  values  are 
assigned  to  the  voxels  whose  intensities  are  far  from  the  centroid. 


To  solve  (I),  we  take  the  first  derivatives  of  J  with  respect  to  Uki,v k  and  set  them  equal  to  zero.  We  thus 
obtain  two  necessary  conditions  for  J  to  be  at  a  minimum: 
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In  implementation,  matrix  V  is  randomly  initialized, then  U  and  V  are  obtained  through  an  iterative  process 
using  (2)  and  (3).  The  stopping  criteria  of  the  iteration  is  that  the  (Euclidean)  distance  between  the  current 


centroids  and  the  centroids  in  the  previous  iteration  is  less  than  a  user  specified  small  number, i.e. ,  \\Vnew  —  V0id ||  < 
e.  After  we  get  V,i.e.,  c  temporal  patterns,  we  select  the  temporal  pattern  with  maximum  initial  enhancement 
as  the  representative  curve  of  the  lesion.  The  corresponding  membership  map,  after  thresholded,  is  the  identified 
region  of  fast  enhancement. 

2.4.  Evaluation 

We  compared  the  signal-intensity  curve  obtained  from  FCM  with  that  from  averaging  the  signal  intensities  over 
the  entire  lesion.  We  extracted  four  features  from  each  curve,  which  was  represented  by  6  time  points,  i.e., 
{so,  si, ...,  S5},  where  Si(i  =  0, 1, ...,  5)  is  the  signal  intensity  at  time  point  i.  Denoting  s*  as  the  maximum  of  the 
six  signal  intensities  and  p  the  time  point  of  s*,the  four  features  are  defined  as: 

maximum  uptake  =  (s*  —  so)/so 
peak  location  =  p 
uptake  rate  =  ( s *  —  so)/(so  x  p) 

{s'-Ss  if  ~  ./  c 

s0x  (5-p) 

0  ifp  =  5 

We  used  linear  discriminant  analysis  (LDA)12  to  combine  the  four  features  in  the  task  of  distinguishing 
between  malignant  and  benign  lesions.  The  differentiation  performance  of  the  LDA  output  in  leave-one-out  cross 
evaluation  was  assessed  using  receiver  operating  characteristic  (ROC)  analysis.13,14 

3.  RESULTS 

Figure  1  shows  a  malignant  case  with  a  breast  lesion  displayed  across  6  slices.  The  membership  map  displayed 
was  thresholded  at  0.4,  i.e.,  the  labeled  voxels  have  membership  values  over  0.4  to  the  cluster  represented  by  the 
solid  curve  in  Figure  1.  The  signal  intensity  time  curve  averaged  over  the  entire  lesion  (dashed  curve)  shows  a 
“plateau”  pattern.  The  FCM  generated  curve,  with  high  initial  enhancement,  (solid  curve)  shows  a  “washout” 
pattern  that  is  characteristic  of  malignant  lesions. 

The  four  features  from  signal-intensity  time  curve  obtained  from  FCM  yielded  an  Az  value  of  0.80,  whereas 
the  features  from  the  signal-intensity  time  curve  obtained  over  the  entire  lesion  yielded  an  Az  value  of  0.68 
(p  value  =  0.007). 


4.  SUMMARY 

Currently  radiologists  manually  select  a  region  within  a  lesion  using  parametric  images.  We  explored  an  unsu¬ 
pervised  learning  algorithm,  which  is  data-driven  and  non-parametric,  with  which  to  analyze  the  heterogeneity 
of  Gd-DTPA  uptake  in  breast  lesions  and  to  automatically  identify  the  high-initial-enhancement  regions.  The 
algorithm  has  the  advantage  of  being  automatic,  fast,  and  efficient. 

The  proposed  FCM-based  algorithm  is  useful  in  automatically  identifying  the  temporal  pattern  with  high 
initial  enhancement  in  CE-MRI  of  breast  lesions  and  improves  the  differentiation  performance  of  enhancement 
kinetics. 


ACKNOWLEDGMENTS 

The  research  was  supported  in  part  by  DOD  Breast  Cancer  Research  Program  DAMD17-03-1-0245  and  USPHS 
grant  CA89452.  M.  L.  Giger  and  U.  Bick  are  shareholders  in  R2  Technology,  Sunnyvale,  CA.  It  is  the  policy 
of  the  University  of  Chicago  that  investigators  disclose  publicly  actual  or  potential  significant  financial  interests 
that  may  appear  to  be  affected  by  the  research  activities. 


1 


2 


3 


4 


5 


Time  (min) 


Figure  1.  A  malignant  case  with  a  3-D  lesion  displayed  across  6  slices.  The  membership  map  displayed  was  thresholded 
at  0.4,  i.e. ,  the  color-labeled  voxels  have  membership  values  over  0.4  to  the  cluster  belonging  to  the  solid  curve.  The 
signal  intensity  time  curve  averaged  over  the  entire  lesion  (dashed  curve)  shows  a  “plateau”  pattern.  The  FCM  curve, 
with  high  initial  enhancement,  (solid  curve)  shows  a  “washout”  pattern  that  is  characteristic  of  malignant  lesions. 
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Abstract — Magnetic  resonance  images  are  often  corrupted  by 
intensity  inhomogeneity,  which  manifests  itself  as  slow  intensity 
variations  of  the  same  tissue  over  the  image  domain.  Such 
shading  artifacts  must  be  corrected  before  doing  computerized 
analysis  such  as  intensity-based  segmentation  and  quantitative 
analysis.  In  this  paper,  we  present  a  fuzzy  c-means  (FCM) 
based  algorithm  that  simultaneously  estimates  the  shading  effect 
while  segmenting  the  image.  A  multiplier  field  term  that  models 
the  intensity  variation  is  incorporated  into  the  FCM  objective 
function  which  is  minimized  iteratively.  In  each  iteration,  the 
bias  field  is  estimated  based  on  the  current  tissue  class  centroids 
and  the  membership  values  of  the  voxels  and  then  smoothed 
by  an  iterative  low-pass  filter.  The  efficacy  of  the  algorithm  is 
demonstrated  on  clinical  breast  MR  images. 

I.  Introduction 

Magnetic  resonance  imaging  (MRI)  has  many  advantages 
over  other  diagnostic  imaging  modalities,  such  as  high  contrast 
between  soft  tissues,  high  spatial  resolution  and  inherent  3D 
nature,  thus  has  gained  wide  clinical  applications.  Breast  MRI, 
for  example,  has  been  widely  investigated  in  the  past  decade  in 
the  detection  and  diagnosis  of  breast  cancer  as  a  complemen¬ 
tary  modality  to  X-ray  mammography[l],  in  assessment  of  the 
localization  and  extent  of  breast  lesions[2],  and  in  monitoring 
tumor  response  to  therapy[3].  Furthermore,  breast  MRI  can 
be  used  for  quantitative  assessment  of  fibroglandular  tissue 
percentage[4]  which  is  a  predictor  of  breast  cancer  risk.  While 
qualitative  visual  assessment  is  often  used  by  radiologists  in 
a  clinical  environment,  computerized  quantitative  analysis  is 
increasingly  needed  to  aid  the  radiologists  increase  both  the 
accuracy  and  the  efficiency  of  the  diagnosis. 

Magnetic  resonance  images  are  often  corrupted  by  inten¬ 
sity  inhomogeneity,  which  manifests  itself  as  slow  intensity 
variations  of  the  same  tissue  over  the  image  domain[5]. 
Such  shading  artifacts  is  the  major  source  of  difficulty  for 
computerized  analysis  such  as  intensity-based  segmentation 
and  quantitative  analysis.  A  number  of  algorithms  have  been 
proposed  for  the  correction  of  spatial  intensity  inhomogeneity. 
Wells  et  al.[ 6],  Guillemaud  and  Brady[7]  developed  a  statis¬ 
tical  approach  based  on  the  expectation-maximization  (EM) 
algorithm  that  simultaneously  estimates  the  bias  field  and 
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segments  the  image  into  different  tissue  classes.  Their  methods 
yielded  impressive  results  on  brain  MR  images  but  has  the 
disadvantage  of  being  computationally  intensive  and  requiring 
prior  knowledge  on  intensity  distributions  of  different  tissue 
classes.  Pham  and  Prince[8]  proposed  an  adaptive  fuzzy  c- 
means  (AFCM)  algorithm  for  fuzzy  segmentation  of  images 
while  compensating  for  intensity  inhomogeneities.  AFCM  is 
robust  in  convergence  because  the  objective  function  to  be 
minimized  has  regularization  terms  that  ensure  the  estimated 
bias  field  is  smooth  and  slowly  varying.  Ahmed  et  al. [9] 
proposed  a  bias-correction  fuzzy  c-means  (BCFCM)  algorithm 
in  which  they  incorporated  a  neighborhood  regularizer  into 
the  FCM  objective  function  to  allow  labeling  of  a  voxel  to 
be  influenced  by  the  labels  in  its  immediate  neighborhood, 
making  the  algorithm  insensitive  to  salt  and  peppernoise.  How¬ 
ever,  they  failed  to  address  that  the  algorithm  may  converge  to 
unwanted  results  without  any  constraint  on  the  bias  field.  Li  et 
a/.  [10]  combined  the  AFCM  and  the  neighborhood  regularizer 
in  BCFCM  and  obtained  promising  results  in  images  with 
high  noise  level.  The  AFCM  based  methods,  however,  are 
computationally  intensive  and  faster  algorithm  is  needed  in 
real-time  clinical  applications. 

In  this  paper,  we  present  a  fast  fuzzy  c-means  (FCM)  based 
algorithm  that  simultaneously  estimates  the  bias  field  while 
segmenting  the  image.  A  multiplier  field  term  that  models 
the  intensity  variation  is  incorporated  into  the  FCM  objective 
function  which  is  minimized  iteratively.  In  each  iteration, 
the  bias  field  is  estimated  based  on  the  current  tissue  class 
centroids  and  the  membership  values  of  the  voxels  and  then 
smoothed  by  an  iterative  low-pass  filter.  The  efficacy  of  the 
algorithm  is  demonstrated  on  clinical  breast  MR  images. 

II.  Methods 

A.  Bias  field  model 

The  observed  MRI  signal  intensity  is  modeled  as  the  “true” 
signal  intensity  multiplied  by  a  spatially-varying  factor  called 
gain  field,  namely, 

Yi=XiGi  V  i€{l,2,...,N} 


(1) 


where  Yit  Xi,  and  Gi  are  the  observed  intensity,  true  intensity, 
and  gain  field  at  the  tth  voxel,  respectively.  N  is  the  total 
number  of  voxels  in  the  MR  image.  The  artifact  can  be 
modeled  as  an  additive  bias  field  by  applying  a  logarithmic 
transformation  to  both  sides  of  (1)[6] 

Vi^Xi  +  pi  V  i  €  {1,2,...,  IV}  (2) 

where  yi,  Xi  are  the  observed  and  true  log-transformed  intensi¬ 
ties  at  the  ith  voxel,  respectively,  and  0i  is  the  bias  field  at  the 
tth  voxel.  By  incorporating  the  bias  field  model  into  a  fuzzy 
c-mean  framework,  we  will  be  able  to  iteratively  estimate  both 
the  true  intensity  and  the  bias  field  from  the  observed  intensity. 

B.  Algorithm 

In  the  image  segmentation  context,  the  standard  FCM 
algorithmfl  1]  is  an  optimization  problem  for  partitioning  an 
image  of  N  voxels,  X  =  into  c  (tissue)  classes 

c  N 

min  {J(U,V-,X)  =  Y,YluPki\\xi-vk\\2}  (3) 

’  k=  1  t  =  l 

subject  to: 

c 

X) Uki  =  i  Vi 

k= 1 

0  <  Uki  <  1  Vfc,  i 

where  U  is  the  partition  matrix  whose  element  Uki  is  the 
membership  of  the  ith  voxel  for  kth  class.  V  is  the  centroid 
vector  whose  element  Vk  is  the  centroid  (or  prototype)  of 
kth  class.  The  parameter  p,  called  fuzzy  index,  is  a  weight¬ 
ing  exponent  on  each  fuzzy  membership  and  determines  the 
amount  of  “fuziness”  of  the  resulting  partition.  The  norm 
operator  ||  •  ||  represents  the  standard  Euclidean  distance.  The 
objective  function  J  is  minimized  when  high  membership 
values  are  assigned  to  the  pixels  whose  intensities  are  close  to 
the  centroid  of  its  particular  class,  and  low  membership  values 
are  assigned  to  the  voxels  whose  intensities  are  far  from  the 
centroid. 

To  incorporate  the  bias  field  into  the  FCM  framework,  we 
substitute  (2)  into  (3).  Then  the  fuzzy  segmentation  with  the 
presence  of  bias  field  becomes  a  constrained  optimization 
problem 

c  N 

min  {Jb{U,V,B-,Y)  =  Y,Y,uUyi-^-v^  <4> 
U’V’B  k= i  «=i 

subject  to: 

C 

X>i  =  l  v* 

fc=i 

0  <  <  1  VA;,i 

where  Y  =  {'tnjfLi  is  the  observed  image,  B  —  {0i}^=i  is 
the  bias  field  image. 


To  solve  (4),  we  take  the  first  derivatives  of  with  respect 
to  Uki,Vk,  and  0k  and  setting  them  equal  to  zero.  We  thus 
obtain  three  necessary  conditions  for  ,Jb  to  be  at  a  minimum. 


*  1 

(5) 

£, I,  <<(!»-« 

Ef.,< 

(6) 

to*  _ 

Pi  -  Vi  ,.p 

2^fc= 1  uki 

(7) 

From  the  first  eye,  an  iterative  scheme  for  minimizing  the 
objective  function  Jb  is  straightforward  by  performing  Picard 
iteration  through  the  above  three  necessary  conditions  for  (4) 
to  be  minimized.  This  is  not  the  whole  story,  however,  because 
B  obtained  from  (7)  is  a  “residual”  image  but  not  necessarily 
be  the  bias  field  image.  A  residual  image  could  always  be 
found  that  would  set  J5  to  zero.  Pham  and  Prince’s  AFCM 
algorithm  solved  the  problem  by  introducing  regularization 
terms  into  the  objective  function  that  ensure  the  resulted  bias 
field  image  being  smooth.  The  regularization  terms,  however, 
make  the  estimation  of  the  bias  field  a  computationally  in¬ 
tensive  process.  Another  solution  is  that  we  estimate  the  bias 
field  by  filtering  the  residual  image  B  in  (7)  using  an  iterative 
low-pass  spatial  filter.  This  filtering  strategy  is  based  on  the 
fact  that  the  bias  field  is  of  low  spatial  frequency  and  the 
assumption  that  other  components  in  the  residual  image  is 
of  higher  frequency.  The  steps  for  our  algorithm  can  then  be 
described  as  the  following: 

1)  Initialize  class  centroid  values,  {ffc}£=1.  Initialize 
{A}ili  with  zeros. 

2)  Update  partition  matrix  U  using  (5). 

3)  Update  class  centroids  V  using  (6). 

4)  Estimate  residual  image  using  (7). 

5)  Filter  the  residual  image  using  an  iterative  low-pass 
filter. 

6)  Go  to  Step  2  unless  the  following  termination  criterion 
is  satisfied: 

H^eio  —  ^o(d||  <  £  (8) 

where  e  is  a  user-chosen  threshold. 

III.  Results 

In  this  section  we  demonstrate  the  efficacy  of  our  algorithm 
by  applying  it  to  clinical  breast  MR  images.  The  images 
were  obtained  using  a  General  Electric  Signa  1 .5-Tesla  clinical 
MR  scanner.  We  set  the  parameter  fuzzy  index  p  =  2,  the 
termination  criterion  e  =  0.001.  The  images  were  thresholded 
before  analysis  so  that  only  tissues  of  interest  were  included 
in  the  computation  and  the  background  was  excluded.  We  set 
the  number  of  classes  c  =  2. 

Fig.  1(a)  shows  a  clinical  MR  image  corrupted  by  intensity 
inhomogeneity.  Fig.  1(b)  shows  the  estimated  bias  field  using 
our  algorithm  and  Fig.  1(c)  shows  the  corrected  image  using 
the  bias  field  in  (b).  The  membership  map  for  dense  tissue 
class  is  shown  in  Fig.  1(d).  Applying  standard  FCM  algorithm 


Fig.  1.  Application  of  the  proposed  algorithm  to  a  clinical  breast  MR 
image:  (a)The  original  image,  (b)  Estimated  bias  field,  (c)  Corrected  image, 
(d)  The  membership  map  for  the  dense-tissue  class,  (e)  For  comparison,  the 
membership  mapping  from  the  standard  FCM  algorithm. 


Fig.  2.  Clinical  breast  MRI  example:  (a)The  orginal  image,  (b)  Corrected 
image.  (c)The  intensity  profiles  of  a  fat  area  as  marked  by  a  line  in  (a),  for 
both  the  original  and  the  bias  corrected  images. 


to  original  image  (a)  without  considering  the  bias  field  effect 
yielded  a  membership  map  shown  in  Fig.  1(e).  Apparently, 
intensity-based  segmentation  could  not  correctly  segment  the 
image  without  intensity  inhomogeneity  correction.  By  visual 
evaluation,  our  algorithm  correctly  estimated  the  bias  field 
and  dramatically  improved  the  image  quality  and  segmentation 
accuracy. 

Fig.2  presents  another  breast  MRI  example.  Fig.2(a)  and 
(b)  are  the  original  and  bias  corrected  images,  respectively. 
We  selected  an  area  mainly  consisting  of  fat  tissue  and 
plotted  the  intensity  profiles  for  both  the  original  and  the  bias 
corrected  images,  as  shown  in  Fig.2(c).  From  the  intensity 
profile  of  the  original  image,  we  can  see  that  the  intensity 
inhomogeneity  in  our  clinical  database  can  be  as  large  as 
40%.  From  the  corrected  intensity  profile,  the  algorithm  has 
successfully  removed  the  bias  field. 

In  both  examples  shown  above,  the  algorithm  converged 


within  1-2  seconds  on  a  PC  AMD  Athlon  with  1.2  GHz  CPU 
speed. 

IV.  Discussion 

We  are  developing  a  computer-aided  diagnosis  (CAD)  sys¬ 
tem  for  breast  MR  imaging.  Bias  field  correction  is  a  necessary 
preprocessing  step  for  subsequent  computerized  quantitative 
analysis.  And  since  the  CAD  system  will  ultimately  be  used  in 
a  clinical  environment, it  must  run  efficiently.  Due  to  the  high 
spatial  resolution  and  3-D  nature  of  MR  images,  we  are  also 
using  MR  images  to  estimate  the  percentage  of  dense  tissues  in 
the  whole  breast,  which  is  an  indicator  for  breast  cancer  risk. 
The  risk  assessment  application  requires  segmentation  of  the 
images  into  different  tissue  classes  in  the  presence  of  intensity 
homogeneity.  All  these  applications  motivated  us  to  develop 
a  reliable,  fast,  and  robust  algorithm  to  solve  the  bias  field 
correction  problem. 

Our  preliminary  experiences  with  the  proposed  algorithm 
showed  that  it  is  a  promising  method  for  intensity  inhomo¬ 
geneity  correction  and  fuzzy  segmentation  of  MR  images. 
Our  work  in  progress  includes  optimization  of  the  current 
implementation  and  evaluation  of  the  method  with  more  data. 
The  noise  sensitivity  of  the  algorithm  will  also  be  investigated 
and  the  neighborhood  regularizer  as  proposed  by  Ahmed[9] 
might  be  incorporated  into  the  current  framework  to  improve 
the  segmentation  accuracy  on  noisy  images.  Finally,  the  cur¬ 
rent  version  of  algorithm  works  for  2-D  images  and  it  is 
straightforward  to  generalize  to  3-D  images. 

V.  Conclusion 

We  have  presented  a  fast  algorithm  based  on  FCM  for 
intensity  inhomogeneity  correction  and  segmentation  of  MR 
images.  The  algorithm  was  formulated  by  introducing  the 
bias  field  model  into  the  FCM  objective  function  which  is 
then  minimized  iteratively.  In  each  iteration  step,  the  bias 
field  was  estimated  based  on  the  current  tissue  class  centroids 
and  the  membership  values  of  the  voxels  and  then  smoothed 
by  an  iterative  mean  filter.  The  efficacy  of  the  algorithm  is 
demonstrated  on  clinical  breast  MR  images. 
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1.  Introduction 

The  purpose  of  this  study  is  to  investigate  the  use  of  a  fuzzy  c-means  (FCM) 
algorithm  for  the  assessment  of  3-D  tumor  extent  from  contrast-enhanced  magnetic 
resonance  images  (CE-MRI)  of  the  breast.  Our  database  contained  121  cases  (77 
malignant  and  44  benign  cases). 

2.  Methods 

For  each  case,  we  initially  delineated  a  cubic  volume  of  interest  (VOI)  that  included 
the  consecutive  slices  about  the  3-D  mass  lesion.  For  each  voxel  in  the  ROI,  the  signal 
intensities  at  six  time-points  were  automatically  divided  by  the  precontrast  signal 
intensity,  yielding  enhancement  curves  .  The  standard  fuzzy  c-mean  (FCM)  clustering 
algorithm  was  then  employed  to  classify  each  voxel  in  the  ROI  into  two  classes:  lesion 
and  non-lesion.  FCM  outputs  two  enhancement  prototypes  associated  with  the  two 
classes  and  the  fuzzy  memberships  of  each  voxel  to  these  two  classes.  The 
enhancement  prototype  that  has  higher  postcontrast  enhancement  represents  the  lesion 
class.  The  voxels  with  larger  membership  values  to  the  lesion  class  were  labeled  as  1 
and  all  others  were  labeled  by  0.  A  3-D  region  labeling  was  then  performed  on  the 
binary  ROI  and  the  largest  region  was  taken  as  the  segmented  lesion.  Finally,  a  3-D 
hole-filling  operation  was  performed  to  include  necrotic  tissue  that  may  have  been 
excluded  from  the  lesion  in  previous  steps  due  to  low  enhancement.  The  results  were 
compared  with  those  of  an  experienced  radiologist  using  the  overlap  measure.  Overlap 
is  defined  as  the  volume  of  the  intersection  of  the  radiologist-delineated  lesion  and  the 
computer-segmented  lesion  divided  by  the  volume  of  their  union. 

3.  Results 

The  mean  overlap  in  3-D  between  the  results  of  the  computerized  segmentation  and 
the  results  of  the  radiologist  outlining  was  found  to  be  0.60  with  a  standard  deviation  of 
0.14.  At  an  overlap  threshold  of  0.4,  91%  of  the  121  cases  were  correctly  segmented. 
The  FCM  method  performed  better  than  our  prior  volume  growing  method  which 
yielded  84%  correct  segmentation  at  the  overlap  threshold  of  0.4. 

4.  Conclusion 

The  proposed  method  is  promising  for  accurate,  consistent,  and  efficient  assessment  of 
breast  tumor  extent  in  contrast-enhanced  MR  images. 


